﻿using System;
using static System.Math;

namespace GeoFits
{
    /// <summary>
    /// 几何图形类：平面
    /// </summary>
    public class Plane : Geometry
    {
        #region 变量和属性
        private Point _pt;
        /// <summary>
        /// 设置或获取平面上的一个点
        /// </summary>
        /// <value>设置或获取平面上的一个点</value>
        public Point PointOnPlane
        {
            get { return _pt; }
            set { _pt = value; }
        }
        private Vector _nor;
        /// <summary>
        /// 设置或获取平面的法向量
        /// </summary>
        /// <value>设置或获取平面的法向量</value>
        public Vector Normal
        {
            get { return _nor; }
            set
            {
                if (value.Length == 0.0)
                {
                    throw new VectorZeroLengthException();
                }
                _nor = value.Normalize();
            }
        }
        /// <summary>
        /// 直线方程系数a
        /// </summary>
        /// <value>直线方程a*x+b*y+c*z+d=0中的系数a</value>
        public double a
        {
            get { return Normal.X; }
        }
        /// <summary>
        /// 直线方程系数b
        /// </summary>
        /// <value>直线方程a*x+b*y+c*z+d=0中的系数b</value>
        public double b
        {
            get { return Normal.Y; }
        }
        /// <summary>
        /// 直线方程系数c
        /// </summary>
        /// <value>直线方程a*x+b*y+c*z+d=0中的系数c</value>
        public double c
        {
            get { return Normal.Z; }
        }
        /// <summary>
        /// 直线方程系数d
        /// </summary>
        /// <value>直线方程a*x+b*y+c*z+d=0中的系数d</value>
        public double d
        {
            get { return -(a * _pt.X + b * _pt.Y + c * _pt.Z); }
        }
        #endregion

        #region 构造函数
        /// <summary>
        /// 构造一个默认平面
        /// </summary>
        /// <remarks>构造一个默认平面，默认平面为XOY平面。</remarks>
        public Plane()
        {
            PointOnPlane = new Point();
            Normal = new Vector(0, 0, 1);
        }
        /// <summary>
        /// 用平面上的一个点和平面的法向量来构造平面
        /// </summary>
        /// <param name="pt">平面上的一个点</param>
        /// <param name="v">平面的法向量</param>
        /// <exception cref="System.Exception">当法向量为0向量时抛出异常。</exception>
        public Plane(Point pt, Vector v)
        {
            PointOnPlane = pt;
            if (v.Length == 0.0)
            {
                throw new VectorZeroLengthException();
            }
            Normal = v;
        }
        /// <summary>
        /// 用平面方程a*x+b*y+c*z+d=0的四个参数来构造平面
        /// </summary>
        /// <param name="a">a*x+b*y+c*z+d=0中的a</param>
        /// <param name="b">a*x+b*y+c*z+d=0中的b</param>
        /// <param name="c">a*x+b*y+c*z+d=0中的c</param>
        /// <param name="d">a*x+b*y+c*z+d=0中的d</param>
        public Plane(double a, double b, double c, double d)
        {
            Normal = new Vector(a, b, c);
            double absa = Abs(a);
            double absb = Abs(b);
            double absc = Abs(c);
            if (absa >= absb && absa >= absc)
            {
                PointOnPlane = new Point(-d / a, 0, 0);
            }
            if (absb >= absa && absb >= absc)
            {
                PointOnPlane = new Point(0, -b / d, 0);
            }
            if (absc >= absa && absc >= absb)
            {
                PointOnPlane = new Point(0, 0, -d / c);
            }
        }
        #endregion

        /// <summary>
        /// 直线与平面的交点 
        /// </summary>
        /// <param name="lp">直线的参数</param>
        /// <returns>两者的交点</returns>
        /// <remarks>注意：如果两者没有交点或者有无穷多个交点，那么将出现异常。</remarks>
        /// <exception cref="System.Exception">
        /// <para>当直线与平面平行时，抛出异常：平面与直线无交点！</para>
        /// <para>当直线在平面内时，抛出异常：直线在平面内！</para></exception>
        public Point InsectionWith(Line lp)
        {
            Matrix paraMatrix = new Matrix(3, 3);
            Matrix ansMatrix = new Matrix(3, 1);
            paraMatrix[0, 0] = a;
            paraMatrix[0, 1] = b;
            paraMatrix[0, 2] = c;
            ansMatrix[0, 0] = -d;
            paraMatrix[1, 0] = lp.Direction.Y;
            paraMatrix[1, 1] = -lp.Direction.X;
            paraMatrix[1, 2] = 0;
            ansMatrix[1, 0] = lp.Direction.Y * lp.PointOnLine.X - lp.Direction.X * lp.PointOnLine.Y;
            paraMatrix[2, 0] = lp.Direction.Z;
            paraMatrix[2, 1] = 0;
            paraMatrix[2, 2] = -lp.Direction.X;
            ansMatrix[2, 0] = lp.PointOnLine.X * lp.Direction.Z - lp.Direction.X * lp.PointOnLine.Z;

            Matrix bigMatrix = new Matrix(3, 4);
            for (int i = 0; i < 3; i++)
            {
                bigMatrix[i, 0] = paraMatrix[i, 0];
                bigMatrix[i, 1] = paraMatrix[i, 1];
                bigMatrix[i, 2] = paraMatrix[i, 2];
                bigMatrix[i, 3] = ansMatrix[i, 0];
            }
            if (Matrix.Rank(bigMatrix) > Matrix.Rank(paraMatrix))
            {
                throw new Exception("平面与直线无交点！");
            }
            else if (Matrix.Rank(bigMatrix) == Matrix.Rank(paraMatrix) && Matrix.Rank(paraMatrix) < 3)
            {
                throw new Exception("直线在平面内！");
            }
            Matrix retMatrix = Matrix.Inverse(paraMatrix) * ansMatrix;
            return new Point(retMatrix[0, 0], retMatrix[1, 0], retMatrix[2, 0]);
        }
        /// <summary>
        /// 点到平面的距离
        /// </summary>
        /// <param name="pt">点</param>
        /// <returns>距离</returns>
        public double DistanceTo(Point pt)
        {
            Vector vt = pt - PointOnPlane;
            double cosAngle = vt.DotProduct(Normal) / vt.Length;
            return vt.Length * cosAngle;
        }

        /// <summary>
        /// 获取平面上离某点最近的点
        /// </summary>
        /// <param name="pt">点</param>
        /// <returns>在平面上与输入点最近的点</returns>
        public Point getNearestPoint(Point pt)
        {
            double dist = DistanceTo(pt);
            return pt + (-dist * Normal);
        }
    }
}
